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ABSTRACT 


The Coppens-Sanders theory for the one-dimensional, nonlinear, 
acoustic wave equation with dissipative term describing the viscous 
and thermal energy losses encountered in a rigid-walled, closed tube 
of large length-to-diameter ratio was applied to finite-amplitude 
standing waves by the use of the Fast Fourier Transform. Camputer 
programs were written to determine the amplitudes and phases of the 
first 255 harmonics. Curves of harmonic distortion as a function of 
the strength parameter were found to be in excellent agreement with 
available experimental data, to agree with the Coppens-Sanders 
perturbation analysis, and to extend the theoretically describable 


régime closer to strengths leading to the formation of the shock front. 
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I. INIRODUCTION 


Theoretical investigations of finite-amplitude waves have been 
characterized by various approximation and perturbation techniques, 
but more significantly by the postulation of various models of the 
absorptive mechanisns [1,2,3,4]. Coppens and Sanders [5] have in- 
vestigated finite-amplitude standing waves for rigid-walled tubes 
of large length-to-diameter ratio by a perturbation method and have 
proposed that the dominant energy losses are produced by the viscous 
effects of the fluid at the walls and the thermal losses due to the 
wall material. They have shown that an approximate wave equation for 
the process is 
CO 

fa 
) (AL* OX. 
er | 
where U>= Se and the subscript n represents the no harmonic 


ie 3 ae 
U,= bY, Stl D+ (1.1) 








of the waveform. The “i frequency component of O, is 


al a OS 
Orn emery, SF ier O, Ne (2. .2) 


where oF is the dissipation parameter [3,5] for the no frequency 





camponent defined by os = ANAS The process has two characteristic 
parameters, one related to the strength of the standing wave and the 
other related to its frequency. 

The purpose of this thesis is extend the Coppens-Sanders approach 
by a "Fast Fourier" analysis technique. By an iterative method, the 
amplitudes and phases of the various harmonics are obtained for dif- 
ferent values of the strength and frequency parameters and are compared 


with available experimental and theoretical results. 


ae) 


Il. PoRORY 


In this section a Fourier-synthesis approach is presented for 
the one-dimensional, nonlinear, dissipative, acoustic wave equation 
based on the theory of Coppens and Sanders as applied to finite- 
amplitude standing waves in a rigid-walled, closed tube of large 
length-to-diameter ratio. 

We begin with Eg. 1.1 which describes the process in question. 


Notice that 


- mee aye 
(ie Pn ~ 7 Polo AY (2.1) 
N= | 


so that Eq. 1.1 may be written as 

















OO 2. 
< es ae 2. 
2 RF Din/ One Pole Ox # | od, 


Assume that p may be written as a Fourier series 





(225) 


€P =) P cost RIL- % Sin nwit+ D, 


where P =Moc on : 
n Onou 


Substitution for a¢ and ©, into the left-hand side of Eq. 2.2 


results in 
(@ 


Ser Ca a ese St a 
on 
Upon substitution for 6 and w_, Exp. 2.4 becomes 
OO n n 


= wae D 2 
\ '-&)2 ee OE lp. 


—_—_— | «see 


3 — = 
rep ax Chote nw dé ot (2.8) 
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c Ze 

Notice that 0a =~ (rk) -@n and O4Pn — fro 

Ds fe eae Pin 
oe that Exp. 2.5 may be written as 

a 
W) Om 3 
Db Hea wera Selo. 
A 2W ° 

\Y\= 
aay W =, + Aw where Aw << Ws and for the fundamental Sz Wr = = 7 


Cy ie 
where the approximation is valid for weak absorption. The approximation 


may be made that 


i) Ses ce [ 2Aw 
_ a ° (257) 


Cae W- 








Define (c,/c,)* =l- 6) and neglect second-order terms; Exp. 2.6 may 


then wee written as 


-2) fs io ar 8 = - ee (2.8) 


\avene®. 
ae right-hand side of Eq. 2.2 contains the quantity om This 











involves multiplication of two infinite series. Terms of the form 
cos p(L-x) sin (qut + Naf for p # gq may be neglected [5], so that pe 


may be spproxinated by 
fp = gel y B, cos RiL-x] sit 


2. 2 
Then we have = (nk) and Eq. 2.2 may be written as 
“% 


OD 
DA —|+ —. 0 2s cosn RlL-x x fsirtneats<ple 10) 
N=\ 


_ Mb p, cosnieiL- x|sin|nat +7} 


Y= \ 





o (279) 











ihe 


If attention is limited to the rigid end at x = L, Eq. 2.10 may be 


considerably simplified to 


a | ARES | 
) |b =e Sa finGust +P) — J coufutsd ig = 


= na 
Peis | 
= ZE)_ Bo sin(rat tf) 





This becomes 


“ Aasin (nuit +P, +6), ) (2.12) 





_ Mb . 
valley. H..Basin(nut +) 
where N=\ 
aus 
@ 2) 2 (2.13) 
y= o, } + a a 7 
a ro, 
and 
| 
-| ” V2 
6}. =a ieee ONE 7. aa : (2.14) 
I+ Rt aS” 
Yr P= | 


The parameters relating to the strength and frequency response of the 
standing wave are Mb/6é 1 and ao , respectively. 

An iterative solution to acl may be obtained by successive 
application of Egs. 2.12, 2.13, and 2.14 to an initial set of AL's 
and |S: The necessary procedure is shown as a flow diagram in Fig. 
1. A discussion of the Fast Fourier Transform and the reasons for its 
use herein are found in Appendix A. Notice that any D.C. terms which 


appear in the procedure have to be set equal to zero since the theory 


excludes these terms. 
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A source of canputational difficulty arises when the frequency 





2 ; ; 
parameter 7% becomes unity. For this case Eq. 2.13 reduces 
1g al 


to Hy =\/n/2, a result which causes the computer solution to become 
unstable because of the strong excitation of the higher harmonics. 
The reason for this difficulty may be seen as follows: the resonance 


frequencies for the classical nondissipative case are defined by 


° T e e ® e 
w= nwo, where WF OT: For linearized propagation in a duct, 


resonance frequencies for exciting each harmonic would be given by 


WO ! = eee eal — = fore) 
a nw nm CLs where CLES (1 6/2) As n > ~, we have the 


é 
limits C7 Sy and w' = Woe For the case of large n then, since 


Be ily, + Aw, the value of the response parameter corresponding to the 


; , 2A 
predicted resonances in the linearized theory is given by 7— ; = 1; 


el 
thus, the higher harmonics are optimally excited and problems will 





arise in the computer program, which admits only the first 255 harmonics. 
Notice that the bulk losses have been neglected up to this point. 


If these terms are retained, may be rewritten [5] as 


go .) 8 ~~ | 
== ligase EN tO : (25) 
3 Ls 


In most experimental configurations the value ofe is usually so small 
that it may be neglected when compared to 6/4n: Injection of an 


unrealistically large value of © does alleviate, to some extent, the 





instabilities in the neighborhood of ao =l . Equations 2.13 
1 
and 2.14 became 
| 
a row | re 
oe 
a nl+ 5 + rae +éEN (2.16) 


Jey 


ave. (a 
‘gua We (2zali?) 


= 
ul iz rA\ 
lls As | eg il ae ts) 
and He would have the asymptotic value of l1fe nas n+», The 


effectiveness of this modification will be discussed in more detail 


in Section III. 
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Initialize A ,o, Mb Zhu 
i 1 w 6 





, Correction Factor 


Calculate Hi, 8, 


Calculate Ps = "a5 





-AL Correction Factor 








| = + pes = 0 
i oe (ce ae >) Correction factor 


Punch cards and print output for AL, ¢, 


Non-computer step: Have A_ and >. reached stable values? 





Recycle with punched cards 


Flow Diagram for Solving Eq. 1.1 


Figure 1 
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III. RESULTS, DISCUSSION, AND CONCLUSIONS 


The computer program written to solve Eq. 1.1 using Egs. 2.12, 
2.13, and 2.14 gave results consistent with the theoretical pre- 
dictions of Coppens and Sanders [5] and Ruff [10]. An additional 
perturbation solution that was performed to insure that the program 
was operating properly may be found in Appendix B. An initially 
sinusoidal pressure waveform when used with a correction factor of 


unity gave rapid convergence to a distorted waveform except that in 
Z2Aw 
WO] 

became unstable. This region of instability gradually increased as 





the neighborhood of response parameter the solution 

the strength parameter was increased. For these cases the modifications 
defined by Egs. 2.16 and 2.17 were used with « =0.00l and allowed a few 
more values of P/?y to be computed, but was not successful in gener- 
ating stable results for the higher strengths. The excellent agree- 
ment with Coppens and Sanders is illustrated in Fig. 2 for the strength 
curves. Figures 3, 4, and 5 illustrate both the area where the solution 
deteriorated and also the frequency characteristics of the various 
harmonics. 

Since the only difficulties encountered were in the vicinity of a 
response parameter of unity, it seems plausible that the instability is 
a result of the behavior of the higher harmonics which become important 
in this region. Recall that a response parameter of unity corresponds 
to the optimum excitation of the higher harmonics and that of necessity 
a truncated series has been used. Then the higher harmonics that are 


retained are affected by the absence of the harmonics that have been 


20 


discarded and react by growling in amplitude to absorb the energy that 
would have been dissipated by those terms that have been amitted. The 
next refinement of the theoretical development of this research would 
be to cause the retained harmonics to behave as if the neglected higher 
harmonics were present. 

Since a solution using iterative processes is easily adapted to a 
computer simulation it is suggested that this method be continued. The 
use of a time-sharing system is recommended since it allows close control 
of the program and saves a great amount of time. 

It has been shown that a Fourier-analysis approach to the Coppens- 
Sanders theory for the one-dimensional, nonlinear acoustic wave equation 
is readily adapted to computer simulation and extends the region of 
validity to the pre-shock régime. The series truncation appears to be 
the only obstacle to obtaining a solution valid for the entire spectrum 


of strength and response parameters. 
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APPENDIX A 


The Fast Fourier Transform 


The Fast Fourier Transform [6,7,8,9] provides an effective means 
for estimating the Fourier coefficients of a set of evenly-spaced 
discrete data. The Fast Fourier coefficients differ only slightly 
from the conventional Fourier coefficients. Construction of wave- 
forms is of prime interest for the higher strength parameters, where 
the Fast Fourier Transform has the distinct advantage of eliminating 
the Gaussian overshoot characteristic of the conventional Fourier 
Transform. This property is illustrated in Figs. 7, 8, and 9. The 
waveform of Fig. 7 was sampled to provide 512 evenly-spaced discrete 
data. This was the input to the procedures available on the IBM-360 


Camputer that perform Fourier and Fast Fourier analyses. The coef- 


ficients from the Fourier analysis were used to construct the waveform 


of Fig. 8, while those of the Fast Fourier analysis were used for Fig. 
9. From the standpoint of computer utilization, the Fast Fourier 
Transform has the advantages of requiring less space and less time 


and was used for the analyses in this research. 


2/ 
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Waveform for Analysis 


Figure 7 
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Reconstructed Waveform Using Fourier Transform 


Figure 8 
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Reconstructed Waveform Using Fast Fourier Transform 


Figure 9 
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APPENDIX B 
Perturbation Solution of Eq. 1.1 


Assume as a first approximation that the solution to Eq. l.l 


consists of only the fundamental which may be written as 
eo ‘ 
fp, = PoGo Meos &(.-%) sin wi. ‘ (Bd) 


Then Z 


2 
B + aM | Es cos2 &(L-x) QA QDies' © 








(B.2) 
COS 2 (| —x) cos2ut| . 
so that the right-hand side of Eq. 2.2 is, at x =L, 
2/2 2 See 
b 3@) SBSIONC M Se -i+coszut (B.3) 
Bas ox 


The squaring of the sin(wt) gives a D.C. term and cos(2uwt). 
The D.C. terms may be excluded [5] so that the second harmonic required 


on the left-hand side of Eq. 2.2 may be assumed to be of the form 


62> rs USS 2® (L-y) s\n (2 wit + Po) (B. 4) 


> ) 
Then -(o" +O > )Pe may be written as 


4% a in + 2.) Ry cosze(L-x) Sve Quit+he) 


ees 


+2VN2Z' BS, cos2Rl_- cout +h \e9 


=| —> , SO that at x = L, Exp. B.5 may be written as 
| 





For Aw = 0, GQ) A 
Be. 





(B.6) 
? 


| / \ Ze }sinf2ust -¢,)-¥ V2" Losleust +2,] 





—4P, 8] 


aL 


which when simplified is 


= 


; (Bey) 
Equating Exps. B.3 and B.7 


= P, R° 2 | 168 Sin 2ut+he pistol =i] 





a PR O 768 Si Al2wt tT Do +9A7.5° 


(B.8) 


= bec. MRE +cos 2wt] 


which may be solved for P, and ¢,. For i = Oca 


2 Sy 


Oe 2 Spey M and One 39radians. 


This agrees with the predicted values of Coppens and Sanders [5] 


and Ruff [10] and the results of the computer program written to solve 


Eq. 1.1. 


SZ 


COMPUTER PROGRAM 


This section contains the computer program used to solve Eq. 1.1 


and a list of the symbology pertinent to the program. 


AAMP(N) =A__) 
FEE(N) = ¢ 5 
BAMP(N) = B-1 
GAMMA (N) = P-1 
H (N) = 785 
THETA (N) = 65 
P (N) = amplitudes of p/Mo.c.* for increment N 


R(N) = amplitudes of (p/Mo.c.*)* for increment N 
ARG = wt/5ll 

RESPM = 20 w/') 9 4 

STRPM = Mb/s, 

CORR = Correction factor 


Computer Program Symbology 
Table I 


‘Spike 


WAWIANA 


AVI oe) 


Silane, 


VAD AOS 


Computer Program for Solving Eq. 1.1 


MENSTCN AASMP(257),F°FF (257) ,Be M06 257),62MMA(257), 
LH(227) THETA (257), TNV(ES) 9S (54) 92 (506) 58 1690) 

1=3. 141593 
AEGS2ePT/E11, 


THE TIALI7ZE FESOM,STRPM,COER, AL MOUN) ~FEE(N). 


PERO (Gee 2 SAM yore Cae 
2 FORME T( 3F 6.3) 

AENN(S,42) (APME(T) ,FEE(T) ,f=1,257) 
4 ENRMAT(3(4X%,27F10.7)) 


TASS CTION CALGULATES Tee ere eee 1 RS tC eee 
THE PHS S baegG eerie | bao ones 
MOD. t=29256 
J=T=-1] 
XN= J 
RaTIC=-(1.0/S0R TOXN) ) 
MY Va eer Poe 
SOE TiN eZ el 36 
HE TI=HAL.C/SOR TC DUMMY] ADIIMMY LER AT ICR RATIC) 
EQUATION 2.13. 
THE TE (TSA TAAQ( PATIO NUMMY] ) 


L (CON TENUE 
eS LOO OR NE c= 152 7 


THTS SECTION CALCULATES THE AMPLITUDEtS BAMP(N) AND THE 
PHESES GAMMATA) CORRESPONDING TO THE CONARE OF THE 
PPFESSURE, 


Sel elem Lege lee 
Sth = 0 oC 
AN 12 J=1,257 
S(TP=HSO(T)+FAAMP( JPESIN(AREX( J-1) *(I-LIFFER(S) ) 
ere COR eh) e 
R(T) )=( OCCT) *O(T)) 
Pi Co iNet 
£5 INPUT CONTAINS TRE AMELITUNES, FOR INCREMENT 
wy yan enor ANPING TO THE SQUARE OE THE OR ESSIRE, 
CALL RHARM (R,R,INV,S,TFERG) 
AS MUTSUT, PIT) CENTAINS THE FAST FOURIEF CO-FFIC IES 
m= THE SQUAPE OF THR PSFSSURS: (80)/72, BC=H=%- 4135 Sie 
(AN) /2, BN=O0, 
Te ANY eee e TEP CREFEICIEN S SF THE CRS INE TERMS sane 
Ti ene COERPETCIENTs OF THE SINE 7 ESMs. 
PAY 13 T=3,511,2 
WPTTE FACH CCMENNFNT GF THE SQUARE OF THE ORESSIIRE £S 
R&MO STINO(NC INCREMENT ) + GEMMA(N)) SY WS ING THE FAST 
PRURTER CHEE“ TeIrNIS. 
J=OT+1) 72 
BAMP( JP=SORT(R(O LT 41) *2 (07 +1) 4807) R07) ) 
GAMMA ( JP=H=ATANZOER(T) ~RCT+1)) 
13 CONTINUE 
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8) 
te 
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“> 
Nw J 


TEES SECT PON ADJUSTS THE LEFT =oaaNe 3 


E OF ENUATION 2.12 
TOWARNS THE RIGHT —H SIDE 


am 
IO 
Ny 
— 


A&MP(T)) 
Ze CER 
SED. s TERMS THy7F RG. 


PPTNT END PUNCH CUT PINT, 


I 
I ( 
3 FNP MS T( 
1] CANTINY 
Y COMNTING 
WRITE ( 
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